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We report on our calculation of the nucleon axial charge g\ in QCD with two flavours of dynam- 
ical quarks. A detailed investigation of systematic errors is performed, with a particular focus on 
contributions from excited states to three-point correlation functions. The use of summed opera- 
tor insertions allows for a much better control over such contamination. After performing a chiral 
extrapolation to the physical pion mass, we find g A = 1.223 ± 0.063 (stat)lg Qgg (syst), in good 
agreement with the experimental value. 



I. INTRODUCTION 

Lattice simulations of Quantum Chromodynamics 
(QCD) have, by now, reached a stage which allows for 
first-principles determinations of many hadronic proper- 
ties, with overall uncertainties at the percent level pQ. 
While systematic errors for quantities such as quark 
masses, meson decay constants and form factors appear 
very well controlled, the situation regarding properties 
of the nucleon is less satisfactory. For instance, lattice 
calculations have so far failed in reproducing the well- 
known experimental findings on nucleon structure (see 
[3] [3] for recent reviews) . A prominent example is the 
axial charge, g A , of the nucleon. Lattice results for this 
quantity lie typically 10 — 15 % below the experimental 
value [4HIB]. What is even more worrying is the absence 
of any tendency in the lattice data which would indicate 
that the gap is narrowing as the pion mass is decreased 
— in fact, the opposite trend is often observed. The most 
likely explanation is that systematic effects are not fully 
controlled. What is lacking, therefore, is a benchmark 
calculation of a quantity which describes basic structural 
properties of the nucleon, and for this purpose the axial 
charge is an ideal candidate: (1) it is derived from a ma- 
trix element of a simple fermionic bilinear which contains 
no derivatives, (2) the initial and final states can both 
be considered at rest, and (3) its definition as an iso- 
vector quantity implies that contributions from quark- 
disconnected diagrams are absent. 

In this paper, we report on our results for g A address- 
ing in detail all sources of systematic errors, such as lat- 
tice artefacts, finite-volume effects, and chiral extrapola- 
tions. We specifically focus on the problem of a system- 
atic bias arising from excited state contributions in the 
relevant correlation functions. To this end, we apply the 
method of summed operator insertions, which helps to 
control any such contamination. 



II. SIMULATION DETAILS 

Our simulations are performed with Nf = 2 flavours of 
0(a) improved Wilson fermions and the Wilson plaque- 
tte action. We stress that excited state contamination is 
an important issue for lattice simulations with any num- 
ber of dynamical quarks. Hence, the question whether 
estimates for g A may be biased can be adequately ad- 
dressed in two-flavour QCD. In particular, there is ample 
evidence [T] that there are no discernible differences be- 
tween QCD with Nf = 2 and N{ = 2 + 1 flavours at the 
few-percent level. Therefore, the observed gap between 
previous lattice estimates of the axial charge and its ex- 
perimental value is by far too large to be explained by 
the presence or absence of a dynamical strange quark. 

We use the non-perturbative determination of the im- 
provement coefficient c sw from ref. |17j . Table[l]contains a 
compilation of lattice sizes and other simulation param- 
eters, including the pion and nucleon masses in lattice 
units. All listed ensembles were generated as part of the 
CLS initiative, employing the deflation-accelerated DD- 
HMC algorithm [TBI US] . Quark propagators were com- 
puted using Gaussian-smeared source vectors [20] sup- 
plemented by HYP-smeared links [2TJ . The smearing pa- 
rameters were tuned to maximize plateau lengths for ef- 
fective masses for a variety of channels. On each ensem- 
ble we collected between 150 and 250 highly decorrelated 
configurations. Up to eight sources, equally spaced in 
the temporal direction, were used to reduce statistical 
fluctuations in correlation functions. In this way, we per- 
formed between 280 and 1700 individual measurements 
on our ensembles. The lattice spacings were determined 
using the mass of the il baryon as described in [52] . As 
we are in the process of supplementing the set of our en- 
sembles, estimates of the lattice spacing will be updated 
in the future. Correlation functions were computed us- 
ing the same smeared nucleon interpolating operators at 
the source and sink. For three-point functions we em- 
ployed the improved axial current which is related to its 
renormalized counterpart via [23j 

(A*) = Z A (1 + b A am q ){A^ + ac A d^P), (1) 
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where and P denote the local axial current and pseu- 
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Run L/a 


ft k amn^ 


amN 


m^L 


A^cfg 


iVsrc 


A2 


32 


5.2 0.13565 0.2424( 4) 


0.592(4) 


7.73 


144 


4 


A3 


32 


5.2 0.13580 0.1893( 5) 


0.531(4) 


6.06 


265 


4 


A4 


32 


5.2 0.13590 0.1454( 7) 


0.481(6) 


4.65 


199 


4 


A 5 


32 


5.2 0.13594 0.1249(14) 0.469(8) 


4.00 


212 


8 


E3 


32 


5.3 0.13605 0.2071( 6) 


0.510(3) 


6.63 


139 


2 


E4 


32 


5.3 0.13610 0.1934( 5) 


0.497(3) 


6.19 


162 


8 


E5 


32 


5.3 0.13625 0.1439( 6) 


0.420(3) 


4.60 


168 


8 


F6 


48 


5.3 0.13635 0.1036( 5) 


0.382(5) 


4.97 


199 


4 


F7 


48 


5.3 0.13538 0.0886( 4) 


0.334(8) 


4.25 


250 


4 


N4 


48 


5.5 0.13650 0.1358( 3) 


0.351(2) 


6.52 


150 


4 


N5 


48 


5.5 0.13660 0.1090( 3) 


0.320(3) 


5.23 


150 


4 



TABLE I. Simulation parameters, pion and nucleon masses 
for all ensembles used in this paper. The temporal extent of 
each lattice is twice the spatial length, T = 2L. 7V c f g is the 
number of configurations per ensemble, while iV src denotes the 
number of different sources. 
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FIG. 1. The ratio R{t,t 3 ) at = 5.2 and = 312 MeV for 
several different values of the source-sink separation t s . 



doscalar density, respectively, and m q is the bare sub- 
tracted quark mass. Since g\ was determined from the 
3 rd component of the axial current alone, the contribu- 
tion proportional to d^P vanishes, as the axial charge 
is defined at zero momentum transfer. Non-perturbative 
values for the renormalization factor were taken from 
ref. [24] . while the improvement coefficient 6 a was es- 
timated in tadpole- improved perturbation theory |25j . 
Since the contribution from the improvement term is at 
the sub-percent level in the range of quark masses con- 
sidered, the systematic effect arising from the unknown 
non-perturbative value for 6a will be negligible. 

III. EXCITED STATE CONTAMINATION 

We denote the Euclidean time separation between the 
nucleon source and sink by t s , while t with t < t s marks 
the interval between the source and the axial current. If, 
as in our case, the same smeared interpolating operators 
are applied at the source and sink, the axial charge can 
be determined from a simple ratio, 

R(t,t s ):=°p^, (2) 
L-2(£ s ) 

where C^(t,t s ) denotes the three-point correlation func- 
tion of the local, bare axial current at vanishing momen- 
tum transfer. For large values of t and t s the ratio R(t, t s ) 
yields directly the bare axial charge, i.e. 

fl(t,t s )*' (t8 -4 >0 5 ba r c + 0(e -A t)+0(e -A( tB - t)); (3) 

where A denotes the gap between the masses of the nu- 
cleon and its first excitation. The axial charge is usually 
extracted by fitting R(t 1 t s ) to a constant. Due to the ex- 
ponentially increasing noise-to-signal ratio in correlation 
functions of the nucleon, typical values of t s are of the 



order of 1 fm. To guarantee a reliable determination of 
5a, excited state contributions in eq. fcfy must already be 
sufficiently suppressed for t, (t s — t) <0.5 fm. 

The lowest-lying multi-particle state in the nucleon 
channel consists of one nucleon and two pions forming 
an S-wave. Therefore, assuming that nucleon and pions 
are only weakly interacting, one expects the gap A to 
be proportional to the pion mass, since the mass differ- 
ence to the Roper resonance amounts to about 500 MeV. 
The same argument applies if one nucleon and one pion 
form a P-wave, provided that the non-zero momentum 
induced by the box size is small enough. It is then clear 
that excited states may increasingly distort the results 
for <7a as the physical pion mass is approached. A typi- 
cal situation is shown in Fig.[T] as t s is varied from 0.8 fm 
to 1.26 fm, the ratio R(t,t s ) is shifted by about 10% to 
larger values. Given the rapid degradation of the signal, 
it then remains unclear whether t s ~ 1 fm is sufficient to 
rule out a bias in the result for <?a- The most straightfor- 
ward strategy to address this problem is to include the 
first excitation into the fit ansatzioi: R(t, t s ) (see ref. [2l)] ) 
or to investigate larger values of t and t s [37] . 

Here we present an alternative approach, based on the 
use of summed operator insertions [201 [2U [29] . The key 
observation is that excited state contributions can be 
parametrically reduced when R(t 7 t s ) is summed over t. 
More precisely, the asymptotic behaviour of the summed 
ratio S(t s ) is given by 

S(t s ) :=£jZ(M„) C + t s {. 9A arc + 0(e- At =)}, 

f=0 

(4) 

where the (divergent) constant, c, includes contributions 
from contact terms. By computing S(t s ) for several 
sufficiently large values of t B , the quantity of interest 
can be extracted from the slope of a linear fit. Since 
t s > t, (t s — t) by construction, excited state contribu- 
tions to the slope of S(t s ) are more strongly suppressed 



3 



S(t, 
30 

28 

26 

24 

22 

20 

18 

16 



11 



a = 
- a = 


079 fm ■ » ■ 

nR^fm - 















































































































().() 



0.7 



().« 



0.9 1 
t, [fm] 



1.1 



1.2 



1.3 



FIG. 2. The summed ratio S(t s ) at rrw « 320 MeV for two 
different lattice spacings (ensembles A5 and F6). 



Run 


a [fin] 


m^pVIeV] 




summ 

5a 




A2 


0.079 


603 


2.454(15) 


1.179( 45) 


1.195(28) 


A3 




473 


2.803(23) 


1.256( 52) 


1.256(28) 


A4 




363 


3.309(41) 


1.084(103) 


1.121(42) 


A5 




312 


3.751(77) 


1.382(127) 


1.228(61) 


E3 


0.063 


649 


2.462(12) 


1.212( 49) 


1.195(40) 


E4 




606 


2.561(11) 


1.154( 68) 


1.160(36) 


E5 




451 


2.920(22) 


1.311(105) 


1.184(52) 


F6 




324 


3.683(48) 


1.268( 91) 


1.217(55) 


F7 




277 


3.771(86) 


1.162( 95) 


1.137(37) 


N4 


0.050 


536 


2.581(13) 


1.221( 32) 


1.176(26) 


N5 




430 


2.939(28) 


1.212( 48) 


1.180(37) 



TABLE II. Results for the axial charge determined via 
summed insertions and the conventional plateau method. 
The data are corrected for finite-volume effects estimated in 
HBChPT (see text). 



relative to R(t,t s ). Compared to the standard method 
of computing the latter at a single fixed value of t s , it is 
clear, however, that the approach via summed insertions 
is computationally more demanding. In Fig. [2] we show 
typical fits to the summed ratio S(t s ) which demonstrate 
that the linear behaviour is very well satisfied. 



IV. RESULTS 

We have determined the axial charge by fitting the 
summed correlator S(t s ) to a linear function for 0.7 fm < 
t s < 1.3 fm and multiplying the slope by the relevant 
renormalization factor of the axial current, eq. ([I]). We 
have verified the stability of the method by excluding 
the smallest value of t s from the fit for each ensemble. 
Typically this leads to an increase in the value for <7a, 
albeit with a 1.5 — 2 times larger statistical error. 

In the following we present a detailed comparison 
between the results obtained using summed insertions 
("summation method") with those arising from fitting 
the ratio R(t,t s ) to a constant in t for t s s» 1.1 fm 
( "plateau method" ) . Results are shown in Table[TT| and 
Fig. [3] One observes that the plateau method yields es- 
timates for <7a that mostly lie below the experimental 
value, which is the typical behaviour seen in other cal- 
culations at similar pion masses. Typically, the summa- 
tion method produces results which are higher than those 
from the plateau method, in some cases by up to 10%. 
At the same time, the summation method has larger sta- 
tistical errors. Still, since an increase is observed in seven 
out of eleven cases, while a slightly smaller value was ob- 
tained only for one ensemble, it is unlikely that this can 
be merely attributed to statistical fluctuations. 

In order to investigate the chiral behaviour in detail, 
we have performed chiral extrapolations based on several 
different ansatze commonly used in the literature [7l[5lfT5l 




0.1 0.2 0.3 0.4 0.5 

m 2 n [GeV 2 ] 



FIG. 3. Chiral behaviour of qa extracted from summed inser- 
tions (upper panel) and using the standard plateau method 
(lower panel). Chiral fits of type "A" and "D" (see text) ap- 
plied for m,Tt < 540 MeV are represented by the blue/hatched 
and red bands, respectively. The black point denotes the ex- 
perimental value. 



H3EDI, i.e. 

Fit A: a + /3ml 

FitB: a 1 + f3' ml- ml In ml /A 2 (5) 

FitC: a" + l3"ml-\i'\e- m * L , 

with fit parameters a, f3, a' , . . .. Another ansatz, Fit D, is 
a three-parameter fit, based on the expressions derived in 
Heavy-Baryon Chiral Perturbation Theory (HBChPT) in 
infinite volume |3 1 j. 132] , with three additional low-energy 
constants fixed by phenomenology [14] . Examples are 



4 



Fit A 



Fit B 



Fit C 



Fit D 



plateau 



summation 



1.1 1.2 1.3 9a 



1.1 1.2 



1.:-! 



FIG. 4. Results for gA at the physical pion mass for the 
plateau and summation methods. Solid points refer to a pion 
mass cut at m w < 540 MeV, while open symbols are used to 
denote results from fits across the entire pion mass range. Fits 
A, B and D were applied to the volume-corrected data (see 
text). The vertical lines represent the experimental value. 
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FIG. 5. Estimates for g\ determined from the summation 
and plateacu methods compared to recent results by ETMC 
[H], RBC/UKQCD ([9] for N t = 2, ref. Q2] for N { = 2 + 1) 
and LHPC ]13,. Only statistical errors are shown. The thick 
vertical line represents the experimental result. 



shown in Fig. [3] A simple linear chiral extrapolation 
(Fit A) applied to the data from all three lattice spac- 
ings for which m ff < 540 MeV yields a value for <?a at the 
physical pion mass which agrees well with experiment 
within the statistical uncertainty. A similar statement 
applies to the fit based on HBChPT (Fit D). By con- 
trast, extrapolations of the data determined using the 
plateau method fail to reproduce the experimental value 
by two standard deviations. 

Fit C was introduced in [5] to test whether the widely 
observed underestimates of gA could be a manifestation 
of finite-volume effects. After determining the parame- 
ters a" , f3" and 7", the volume-dependent term propor- 
tional to exp{— m^L} can be subtracted. Indeed, a non- 
zero value for 7" results when fit C is applied to the data 
obtained via the plateau method. A linear chiral ex- 
trapolation, using the fitted coefficients a" and /3", then 
yields an estimate for gA which agrees with experiment 
(see Fig.[4j. However, repeating the procedure for the 
summation method produces a vanishing coefficient 7". 
We conclude that, in this case, there is no need to sub- 
tract any term designed to account for finite-volume ef- 
fects, in order to get agreement with experiment. When 
addressing the influence of excited states it is important 
to realize that such contributions are volume-dependent, 
whenever they are due to multiple-particle states. Thus, 
for a true benchmark calculation of the axial charge one 
must be able to separate finite-volume corrections to <?a 
from volume-dependent excited-state contamination. 

Figure|4] shows a compilation of the chirally extrapo- 
lated gA from the four different fit types. While summed 
insertions invariably produce estimates that are compat- 
ible with experiment, one consistently obtains lower val- 
ues using the plateau method, except for fit C with the 



term containing exp{— m^L} subtracted. These observa- 
tions are stable under variations of the pion mass range, 
as indicated in the figure. 

We now proceed to discussing our final result and the 
estimation of systematic errors. We applied a finite- 
volume correction based on the expression derived in 
HBChPT HI] (see ref. [14] for details on the numerical 
evaluation). Since m n L > 4 and 2fm < L < 3fm the 
resulting shifts are at the sub-percent level for all our en- 
sembles. As our best estimate, we quote the result from 
fit A, applied to the volume-corrected data obtained from 
summed insertions, with a cut of < 540 MeV, i.e. 

g A = 1.223 ±0.063 (stat), (6) 

which agrees with the PDG average [55] of 1.2701(25) 
within the statistical error. By contrast, when the same 
fitting procedure is applied to the results extracted from 
the plateau method, one finds the much lower estimate 
of g A = 1.173 ± 0.038 (stat). 

It is instructive to compare our findings to other re- 
cent results for the axial charge. A compilation is plot- 
ted in Fig. [5] With the exception of the results by 
RBC/UKQCD [TU], our estimate based on the sum- 
mation method is the only one which agrees with the ex- 
perimental value within statistical errors. It it also worth 
mentioning that, in order to achieve agreement with ex- 
periment, a large downward curvature in the data had to 
be separated off in refs. [9, 10 , by applying the procedure 
of fit C. 

We note that, with our current level of statistical ac- 
curacy, no significant dependence on the lattice spacing 
could be detected. For instance, applying fits A-D only 
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to the data at j3 = 5.3 produces a tiny variation, which 
is 10 times smaller than the statistical error. There- 
fore we refrain from quoting a separate systematic uncer- 
tainty relating to cutoff effects. In order to quantify the 
uncertainty associated with the chiral extrapolation, we 
adopted two procedures. First, by applying different cuts 
to the upper limit on the pion mass interval between 470 
and 640 MeV, we observe a variation of ±0.035 relative 
to the central value in eq. ^ . Second, we considered 
the spread among fits A-D as a measure for the uncer- 
tainty relating to the extrapolation, which amounts to 
a downward shift by —0.060. Taking the largest upward 
and downward variations from both methods as the error 
estimate, we arrive at our final result 

<7A - 1-223 ± 0.063 (stat) ±° ° f (syst), (7) 

which agrees with the experimental result at the level of 
6 - 7%. 



butions are parametrically more strongly suppressed rel- 
ative to those encountered in conventional ratios. Our 
findings, summarized in Fig.[4j demonstrate that a much 
better agreement with the experimental value of gx can 
be achieved in this way. On the downside, one must list 
the necessity to compute correlation functions for several 
source-sink separations, as well as the larger statistical er- 
rors associated with the method. However, since excited 
state contamination might be a generic problem for lat- 
tice calculations of structural properties of the nucleon, 
the larger numerical effort seems a worthwhile invest- 
ment. We plan to corroborate our findings by including 
additional ensembles with smaller pion masses and ex- 
tend our studies to other quantities, such as the vector 
and axial vector form factors of the nucleon. For this pur- 
pose, optimised anisotropic smearing functions for non- 
vanishing hadron momenta [341 135] may prove to be a 
useful addition to the technique of summed insertions. 



V. CONCLUSIONS 

The typical source-sink separations in baryonic three- 
point functions can be smaller by up to a factor two 
compared to those used in the mesonic sector. Even for 
t s ~ 1.3 fm it is hard to judge whether or not a significant 
bias due to excited state contributions can be excluded, 
if the standard plateau method is employed without fur- 
ther checks (see Fig. [I]). Summed operator insertions of- 
fer an attractive alternative, since excited state contri- 
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